#!usr/bin/perl -w
use strict;

#run like this: perl searchforTIR.pl -t 4_Nc_vs_5_Nc.blastn.filterd.tab 
#-db Nc_Mariner_flankingseqs.fasta 

use Bio::DB::Fasta;
use Bio::SeqIO;
use Getopt::Long;
#use Bio::Perl;

my $db;
my $table;

GetOptions('d|db:s' => \$db,
			't|table:s' => \$table,
			);
			
			
open(TAB, "<$table")  || die "cannot open file $table: $!\n";
open(DATA, "<$db")  || die "cannot open file $table: $!\n";

my $qstart;
my $qend;
my $seqID; #we use query firstly

while(<TAB>)
{
	next if /^#/;
	my @row = split(/\s+/,$_);
	$qstart = $row[6];
	$qend = $row[7];
	$seqID = $row[0];
}
print "qstart is $qstart and qend is $qend\n";

my @leftend;
my @rightend;
my $rev_com_right;
my @EX_LF;
my @EX_RT;
my @seq;
while(<DATA>)
{
	next if /^>/;
	my $count;
	my @row = split(//,$_);
	
	#this extracts the left flanking region of the sequence
	#print "the $qstart character is:";
	#print "$row[$qstart-1]\n";
	for(my $i=0; $i < 100; $i++)
	{
		push(@leftend, $row[$qstart-1+$i]);
	}
	print "the leftend seq is @leftend\n";
	
	#this extracts the right flanking region of the sequence
	for(my $i=0; $i < 100; $i++) # this stores the reversed sequence of the rightend
	{
		push(@rightend, $row[$qend-1-$i]);
	}
	my $tmp = join("",@rightend) ;
	$rev_com_right = com($tmp); # this function gives the complement sequence
	@rightend = split(//, $rev_com_right);
	print "the reversed-complement rightend seq is @rightend\n";
	
	#this extracts the possible TIR region after the first occurrence of "TA" on the left end
	my $index;
	for(my $i = 0; $i < 100; $i++)
	{
		if($leftend[$i-1] eq "T" && $leftend[$i] eq "A")
		{
			$index = $i+1; 			
			last;
		}
	}
	while($index < 100)
	{
		push (@EX_LF,$leftend[$index]);
		$index++;
	}
	print "the left sequence after first TA is: @EX_LF\n";
	
	
	#this extracts the possible TIR region before the last occurrence of "TA" on right end
	my $index_r;
	for(my $i = 0; $i < 100; $i++)
	{
		if($rightend[$i-1] eq "A" && $rightend[$i] eq "T")
		{
			$index_r = $i+1; 			
			last;
		}
	}
	while($index_r < 100)
	{
		push (@EX_RT,$rightend[$index_r]);
		$index_r++;
	}
	print "the right sequence before last TA is: @EX_RT\n";
	
	#this finds the longest 
	my $longest_com_sub;
	my $length = 0;
	my $l = join("",@EX_LF);
	print "left string is $l\n";
	my $r = join("",@EX_RT);
	print "right string is $l\n";

}



#this function returns the complemented sequence for a sequence 
sub com
{
	my $seq = $_[0];
	$seq  =~ tr/ACGTacgt/TGCAtgca/;
    return $seq;
}
			

